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ABSTRACT 

We present new perturbation theory (PT) predictions in the Spherical Collapse 
(SC) model for the 2-point moments of the large-scale distribution of dark matter den- 
sity in the universe. We assume that these fluctuations grow under gravity from small 
Gaussian initial conditions. These predictions are compared with numerical simula- 
tions and with previous PT results to assess their domain of validity. We find that the 
SC model provides in practice a more accurate description of 2-point moments than 
previous tree-level PT calculations. The agreement with simulations is excellent for a 
wide range of scales (5 — 50/i _1 Mpc) and fluctuations amplitudes (cr 2 ~ 0.02 — 2). 
When normalized to unit variance these results are independent of the cosmological 
parameters and of the initial amplitude of fluctuations. The 2-point moments pro- 
vide a convenient tool to study the statistical properties of gravitational clustering for 
fairly non-linear scales and complicated survey geometries, such as those probing the 
clustering of the Lya forest. In this context, the perturbative SC predictions presented 
here, provide a simple and novel way to test the gravitational instability paradigm. 



1 INTRODUCTION 



The 2-point point correlation function £2 is a well known 
and useful (eg Peebles 1980, 1993) characterization of fluctu- 
ations in the large scale density distribution. It measures the 
covariance of density fluctuations £2(^1, ri) =< <5(ri)5(?*2) > 
and should be translationally invariant and isotropic, i.e, 
£a(n,r2) = £2(»"12) = |»"2 — n|), in an homogeneous and 
isotropic distribution. In fact, for a Gaussian field, £2(^12) is 
all we need to characterize fluctuations, as all higher order 
correlations £j\r(ri, . . . ,Tn) =< S(ri) ■ ■ ■ 5(tn) >c are zero 
(the subscript c stands for the 'connected' or reduced mo- 
ments). Even if the initial conditions were Gaussian, non- 
linear evolution of the density field under gravity intro- 
duces non-Gaussian correlations so that both in the weakly 
non-linear regime (e.g, Fry 1984) or the strongly non-linear 
regime (e.g, Davis & Peebles 1977) we expect £jv — ^ 2 V_1 - 
It is therefore important to characterize and predict in de- 
tail the behavior of higher order correlations. One option 
is to consider £jv itself, but this is complicated because of 
the multimensional space n, . . . , rjv involved. Another op- 
tion is to consider N-order one-point moments or cumulants 
K.N =< S N >, smoothed over a fixed scale R. This has the 
advantage of providing good signal to noise in the measure- 
ments (as we are averaging over multidimensional space) 
but at the cost of losing configuration information. Another 
problem is to take into account the smoothing window when 



we have complicated geometries. For example when the sur- 
vey consists of one dimensional pencil beams, as is the case 
of the density traced by the Lya forest in QSO spectra. This 
is in fact the main physical applications that motivated the 
present study. 

Intermediate between cumulants and N-point correla- 
tions are 2-point cumulants: 

fC pq ( ri ,r 2 ) =< 5(r 1 ) p 5(r 2 ) q >c 



where Kpairi.ra) = IC P 



(1) 

(7-12), in a homogeneous and 
isotropic distribution. These 2-point cumulants have some of 
the advantages and disadvantages of N-point functions and 
N-order cumulants. They do contain information on higher- 
order correlations (K, pq = 0, for p+q > 2 in a Gaussian field), 
but are better suited than cumulants for situations where we 
want to avoid boundary restrictions in complicated geome- 
tries. As mentioned previously, a good example is the Lya 
forest, which traces a narrow one dimensional sample of the 
density distribution along the line-of-sight of the QSO. 

As in the case of the cumulants and N-point functions 
one can use perturbation theory (PT) to predict the values 
of Kpq in gravitational perturbation theory from some given 
(typically Gaussian) initial conditions. Predictions for the 
2-point moments in the PT approach to lowest order were 
first derived by Bernardeau (1996; B96 hereafter). 

This paper is motivated in part by the lack of measure- 
ments of 2-point moments, K, pq , from numerical simulations. 
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As pointed out in B96, numerical results for K. pq are difficult 
to obtain in a regime where they can be directly compared 
with the PT predictions. However, simulation work is an im- 
portant and necessary step, as we need to study both the 
non-linear regime and the range of validity of the PT calcu- 
lations for the 2-point statistics. 

As we will show below this comparison can be most con- 
veniently done by computing the PT predictions for the den- 
sity moments using the Spherical Collapse (SC) model. The 
SC model has been shown to provide a simple and accurate 
way of predicting the 1-point statistics of gravitational clus- 
tering both in the weakly non-linear (Bernardeau 1992, 1994; 
Fosalba & Gaztafiaga 1998, FG98 hereafter; Gaztanaga & 
Fosalba 1998) and fully non-linear (see e.g, Gaztanaga & 
Croft 1999; Gaztafiaga & Scherrer 2001) regimes. 

In this paper we concentrate on the perturbative ap- 
proach. The fully non-linear case and applications to Lya 
forest will be presented elsewhere (Gaztanaga et al. 2001). 
In paper II of this series we present results in redshift space 
and a comparison with clustering in galaxy data. 

This paper is organized as follows. In §0 we introduce 
the 2-point statistics. Perturbative predictions in the SC 
model are given in §H and a detailed comparison to Nbody 
simulations is presented in iM Finally, we summarize our 
results and give our conclusions in §0. 



P°2(Si,S 2 ) « Pf(5i) Pf(<5 2 ) [l + a 6 1 6 2 + O (a 2 )] . (6) 

As we shall see below (see §H) the expansion in Eq.(H) is 
a basic ingredient for generating PT predictions for the 2- 
points moments, as it simply factorizes the 2-point average 
into products of 1-point averages. 

2.2 The Evolved Mass PDF 

2.2.1 Linear Evolution 

Because of gravitational growth, the evolution of the matter 
density field, 8, will change the initial PDF. For small fluc- 
tuations linear theory provides a simple prediction for the 
time evolution: 



S(t,x) = D(t) 6 (x) 



(7) 



where D(t) is the growth factor (equal to the scale factor 
D = a for Q. — 1), and 5q(x) is the initial field. We will de- 
note this linear prediction by Sl- For Gaussian initial condi- 
tions the linear PDF is also Gaussian with a 2-point function 
£2(2), given by scaling the initial value £2(20) by D 2 (t), so 
that: 

£ 2 (t)=D 2 (t)£ 2 (io). (8) 



2 THE 2-POINT PDF AND ITS CUMULANTS 

2.1 The PDF of the Initial Conditions 

In the limit of early times, we assume a nearly homogeneous 
matter distribution with very small fluctuations (or seeds), 
with given statistical properties. We will concentrate on the 
case where the statistics of the initial density are well de- 
scribed by Gaussian initial conditions, which correspond to a 
broad class of physical models to generate initial conditions. 
The 2-point PDF is given by a bivariate normal distri- 
bution for Pg(<5i, fe): 



P&Vuh) = 



exp 



■sE * 



°ii 5 j 



(2) 



2nVdetC 

where C is the inverse of the covariance matrix: 

dj = < Si Sj > i, j — 12. (3) 

In our case we have that en = C22 = o 2 , and C12 = C21 = 
£2(^12) is the 2-point correlation function. Notice that a 2 — 
£2(0). We then have: 

1 



■fl2 



2^ 4 - I 



exp 



a 2 S 2 + a 2 S 2 - 2£ 2 <5i<5 2 



2(<r 4 -£ 



(4) 



As correlations decrease with separation, we typically have 
£2(^12) << o- 2 . In this limit we find: 



p&tfuh) 



P 1 G (5 1 ) P?(5 2 ) 



1 + 



£2(^12) 8182 



(•>) 



where P\(S) is the 1-point PDF for a Gaussian field [j. Thus, 
the explansion parameter is a = £,2(ri2/o~ 4 , and we have: 



2.2.2 Non-linear Evolution 

The (Newtonian) non-linear equations of motion of density 
fluctuations in the matter dominated regime have no exact 
analytic solution (Peebles 1980). However, one can conve- 
niently describe the non-linear evolution of density fluctua- 
tions by making use of the Spherical Collapse approximation 
(SC hereafter). Within this approximation the non-linear 

overdensity p = 1 + 8\\ in Lagrangian space is simply related 
to the linear one by: 



p(q) = G[8L(q)] 



(9) 



where 8l{q) = D8o{q), and we stress that we choose La- 
grangian coordinates q to describe the growth of fluctua- 
tions. The non-linear 2-point probability distribution func- 
tion of matter density fluctuations, P {81,82), is thus related 
to the linear (Gaussian) one by a simple change of variable: 



P( Pl , P 2)=P G (8tSi) 



d8{ 
dpi 


d8k 
dp2 



(10) 



where Sl — <5 _1 [p]. 

As we are interested in the weakly non-linear regime, we 
shall concentrate on the perturbative (or Taylor) expansion 
of the non-linear density fluctuation in terms of the linear 
one: 



p(q) = G[5L(q)]=J2^TVL(q)r 



(11) 



where the v coefficients specify the non-linear coupling be- 
tween density fluctuations in real-space. Note that in the SC, 
these z/s are the equivalent to the angle average (spherically 
symmetric component) of the couplings between modes as 



Note that a Gaussian PDF only makes physical sense, i.e, 8 > 
-1 (p > 0), when the variance is small: a — » 



T As it is more convenient, in the equations below we will use the 
notation 1 + 5 for the non-linear case and 8l for the linear one. 
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given by the PT kernels in Fourier space (see FG98 for de- 
tails). The first such coefficients are (see Bernardeau 1992, 
1994; FG98), 



v-i 



I'A 



1'4 



V5 



£-.- 



682 
189 
446440 
43659 ~ 
8546480 
243243 



3.61 



10.22 



35.13 



(12) 



2.2.3 Lagrangian to Eulerian PDF 

As mentioned before, the above expression Eq.(|10[) for the 
PDF corresponds to the probability distribution of the 
evolved field in Lagrangian space, q. This corresponds to 
a fixed mass, while we want the statistics for a fixed vol- 
ume. The Lagrangian and Eulerian elements can be related 
by imposing mass conservation: 



dSi(x) dS 2 (x) - pi p 2 dSi(q) dS 2 (q), 



(13) 



where p(x) = 1 + 5(x) is the overdensity in Eulerian coordi- 
nates. We then have: 



Pe{Pi,P2) 



P1P2 



Pl(pi,P2) 



(14) 



up to a normalization constant. The Eulerian PDF is there- 
fore 



Pe(pi,P2) 



1 Pg{8£,8Z) 



N 



P1P2 



d6£ 

dpi 


dp2 



(15) 



where N is a normalization constant. 

In the perturbative regime one finds that the PDF in 
Euler and Lagrange space are the same, as the term p\p2 
in the denominator above only introduces changes in higher 
order corrections (i.e, to leading order, p\p2 — 1). The next- 
to-leading order in the expansion of the PDF or its mo- 
ments, the so-called one-loop correction (see Scoccimarro & 
Frieman 1996a, b) has a contribution from the mass conser- 
vation factor, P1P2 — pi, which yields, 



Pe(pi,P2) 



1 P G (St,5^) 



N 



Pi 



dSf 
dpi 


dp2 



(16) 



and which we shall use below for the estimation of the 2- 
point moments of the dark matter density field in the weakly 
non- linear regime Fl. 

2.3 2-point Cumulants 

Consider now a generic field, (f>, which could be either the 
measured flux in a ID quasar spectrum or the mass den- 
sity in 3D space. We will assume that this field has been 
smoothed over some resolution cell A. We define the reduced 
2-point reduced moments or cumulants tC pq of the field <f> at 
two different positions ri and r 2 , by: 



■f This expression apparently breaks the 1^2 symmetry, but 
this is only a trick that will be useful for simplifying the moment 
estimations below (see also the argument in §3.2.2 of B96) 



/Cn 



E< ^(n) p ^(r 2 ) q >c 



d 2 logM[ti,t 2 ] 



dtidh 



(17) 



*l,*2-»0 



where M\ti,t2\ =< exp(<j!>iti) exp(02^2) > is the joint gen- 
erating function of the (un-reduced) moments: 



m pq =< <j)(ri) p '<j>(r 2 ) q >— 



d 2 M[ti,t 2 ] 



dtxdtz 



(18) 



which can also be readily obtained by using the correspond- 
ing 2-point probability distribution: 



m pq =< (j>(ri) p (j)(r 2 ) q > = 



m 



id(f> 2 (19) 



For example, the 2-point function of the field is just 
given by: 



£2(7*12) = lCn(ri 2 ) —< <t)i(t> 2 > c 



(20) 



Of course, as mentioned in the introduction, these 2-point 
cumulants are just a function of the separation 7*12 = r 2 — ri 



between the cells: K, v 



- IC P q(ri 2 ) 



The first reduced moments are: 



A^oo 


- m o 




A^oi = 


= m i 




ICq 2 - 


2 
- a — 


m 2 — 77J 1 


A^03 


- "103 


— 3moim 2 + 2moi 


JC04 - 


= m 4 


— 6m i + 12m i m2o — 3m 20 


ICu - 


= & = 


"In — 771 i 


lCl 2 


= TO12 


— ?Tl01?Tl02 — 27710177111 + 2l71 



(21) 



4m 177103 



and so on. Note that even when we normalise the field so 
that the mean is zero (ttioi = 0), the cumulants are different 
from the central moments in that they have the lower order 
moments subtracted, i.e: 



A^22 = 


= 77122 


- 2t71ii — 7712077102 




A^13 = 


= 77113 


- 37710277111 




K.14 


= 77114 


- 4mo377ln — 67710277112 




K-2A 


= 77123 


- 6mi2T71ii ~ 37710277121 - 


- 7710377120 



(22) 

for 77101 = 0. 

It is interesting to define the following 1-point hierar- 
chical ratios: 

c _ fog 



q> 2 



(23) 



and the corresponding 2-point generalization: 



^"d 



lC pq {r i 2 ) 



ICu KZt q ' 2 ?2( r i2) t 2(p +"- 2) 



p + q > 2 (24) 



At 7"i2 = the 2-point cumulants become 1-point cumulants, 



fC, 



fCi+j and consequently, for the reduced moments, 



c pq — > S p + q . These ratios turn out to be roughly constant 
under gravitational evolution from Gaussian initial condi- 
tions (see Peebles 1980). 



2.4 Discreteness Effects 

For discrete fields, such as fluctuations traced by the galaxy 
distribution or particles in Nbody simulations, we must 
correct the cumulants to account for Poisson fluctuations 
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around the mean density N =< TV > at the smoothing scale. 
Following Gaztanaga & Yokoyama (1993), there is a simple 
relation between the discrete cumulant generating function 
Af D (t) and the corresponding continuous one Af c (t): 



N D [t]=N c 



1], 



(25) 



where N(i) generates the cumulants (or reduced moments) 
of the density field: N(r). The latter is simply related to 
the density fluctuation 8(r), N(r) = 7V[1 + 5(r)]. We can 
trivially extend the above relation to the 2-point case: 



M D [ti,* 2 ]=Af c [ei- 1,4-1] 



(26) 



This relation yields the following corrections for estimating 
the continuous cumulants K,p q from the measured, discrete, 
ones K,„„: 



(27) 



^02 ■ 


yD 
= ^-02 - 


1 






^03 = 


yD 
= ^03 - 


3 y-C 

= >M)2 ~ 

N 


1 


5 


K\\ 


yD 








K-\i - 


yD 
= *-12 ~ 


=/^iij 

TV 






^13 = 


= ^13 _ 


3 „.c 

N 


1 

~W 2 


£n, 


^22 = 


= *- 2 2 - 


2 „c 

= A~12 " 


1 

^2 


£ll! 



TV 



TV 



and so on. In all cases /C pg = K pq (ri,T2) with n 7^ ?*2- Note 
how the two point function K,\\ is not affected by Poisson 
fluctuations, unless the two points are the same r-y = V2 (e.g, 
ICo2). This reflects the fact that Poisson fluctuations are not 
spatially correlated and it only yield contributions when we 
have two or more cells at the same location. 

In the limit when shot-noise dominates the cumulants 
we have that the discrete fields give for Eq.(E3|): 



S 



/c, 



0,9 



(Ko 1 



= 1 



q> 2. 



(28) 



In this same limit Cp°' aBon = 0. For a Gaussian continuous 
field (i.e, neglecting only the clustering of higher orders p + 
q > 2 ), we find: 



Poisson i 

C12 = 1 ; 

Poisson o Poisson A 

C22 = o , c 13 = 4 ; 

Poisson t o Poisson n n 

c 23 = 18 , c 14 = 32 



(29) 



and so on. These values illustrate the fact that Poisson fluc- 
tuations produce artificial non-Gaussianities which mimic 
the hierarchical clustering. We will see below that the hier- 
archical amplitudes that emerge from gravitational growth 
(starting from Gaussian initial conditions) are significantly 
higher for all realistic situations. Nevertheless, the shot-noise 
must be subtracted (using the formulae above) if we want 
an accurate comparison with predictions from gravitational 
instability. 



3 PERTURBATIVE PREDICTIONS 



3.1 2-point Cumulants for the Mass 

We can now use the 2jpoint quasi-linear distribution induced 
by the SC, i.e, Eq.dlq), to estimate m pq and KL pq for mass 
density fluctuations. We first consider the case where the res- 
olution cell A — *• 0. This is a straightforward, rather tedious 
calculation (ideally fitted for algebraic programming, such 
as Mathematica). The first step is to expand the moments 
in a perturbative series. We use Eq.(B) and the expansion 
Eq.dlll) for the non- linear overdensity Q: 



m vq = < Q\ Q\ > = Ipfilqfi + Q Ip,llq,l 
~^~ 9 ° 1/p> 2 ^«> 2 ~~ i^Pfi^lfi + Iq,2lpfi) CT + Ipfilqfi & \ 

+ - a 3 [J p ,3/g,3 — 3 (7 P , 3/9,1 + Iq,3l P ,i) cr 2 + 

9 I P ,iI q ,i a 4 ] + O (a 4 ) ; 

im,» = ({Q-1T 51) (30) 

where we have truncated the expansion to third order in 
ct — ^2(^i2)/c 4 , and the mean < . . . > is taken over the 1- 
point Gaussian PDF Pg{5). Note that in Lagrangian space 
we also need to include a term 1/Q — 1/p to account for 
mass conservation. 



3.1.1 Results from the SC Model 

As discussed in §@, the SC approximation is derived in La- 
grange space, but clustering in the simulation and observa- 
tion samples is measured in Euler space. To normalize ap- 
propriately the moments to Euler space, we shall take into 
account mass conservation, as discussed in §2.2.3. This im- 



plies that the 2-point moments in Lagrange and Euler space 
relate to one another in the following way: 



<PlP2 >E = < Pi P2 >L <l/p> L 



J+K-l 



(31) 



where the sub-indices E and L denote Euler and Lagrange 
space, respectively. The analogous expression for the 1-point 
moments is recovered by setting K = 0. 

Introducing Eq.(Bll) into Eqs.(E2) & (pO[), one can derive 
PT predictions for the 2-point moments to arbitrary order. 
In particular, we find that all the 2-point cumulants JC pq 
(with p, q > 0) will be proportional to the two point function 
£2 = /C11, 



^Cp 9 (»"12) 



(<T 



2Np+q-2 



&(r«)+0(G) 



(32) 



where a 2 — JCn(0) = IC02 — IC20 — °" 2 (A) is the variance at 
the resolution cell. 

To leading order in the two point function, £2 = ICu, 
they are given by: 

AC11 = [l + (3 - 3^2 + v 3 )o- 2 ] £2 (33) 

+ (A -A ^+°(^) 

IC02 = a 2 + (3 - 4^ 2 + i*/ 2 2 + jy 3 ) o 4 + 0{a 6 ) 

K.12 — \2V20- — (2 — 13^2 + 9^2 + 3^3 — 3^2^3 — Va) <T ] ^2 
K.Q3 — ?>U20 
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-2 + 18^2 - —i4 -4i/ 3 



+ l/| + 6i/ 2 *>3 + ? t/ 4 ) a 6 + O (a 8 ) 



Similarly, for the reduced moments (or hierarchical ra- 
tios), Eq.(gJ), we get, 

Cl2 = 2^2 + ( — 2 + V2 + hv\ — 3^3 — v\ — ^2^3 + ^4) O" 2 



+ 



}M 



-^ — V2 + -U-2 + -{—V2 + Vi) 



6 



ti 



CIS 



— 6z/ 2 + 3^3 

+ (6 - 30^2 + 18^2 + 30i/f - 6V2 - 2^3 - 9u 2 u 3 



(34) 



+ 



+ 



%V%Vz 



Qvl 



-v\ — Aua +9^2^4 + -^5 ) O" 2 



6^2 — 3^2 + 24f 2 — 6/^2 



3 2 15 

9^2^3 — 7J^2^3 + -^-V2VA 



—5+6is 2 — 3u 2 — 1v-i — 3u 2 v-s 



£5 



+ 5z^3 + 2^3 + y 



/ 1 2 3 

+ I ^2^3 — X^2^3 — V4, + -V1U4 



$ 



2 + oat) 



(35) 



C22 = 4^2 + ( —4^2 — 6^2 + 24^2 



•4wJ 



6^2^3 — 4i/ 2 ^3 + 4^2^4J CT 



+ 



+ 



-^(4 + ^3) + 6 - 24^2 + lli>f + 22i/| 



4^2 — 6^3 + 6^2 ^3 — 8^2 ^3 — 2^3 



3^4 +6^2 ^4 + 2^5] £2 



2 

4-| - 8^2 + 81/3 + 16i/2 - 6i>| - 14^3 

a 2 



+ 2l^l/3 ~ 2l/ 4 + 8K2I/4I ^ 

+ (2^ - 2l/| - f2 + 2^3 



(36) 



The above expressions, Eqs.(|34|)-(|3q), are the main analytic 
results of this paper. Notice that, to leading order (tree level 
in PT) for Gaussian initial conditions, the following property 
holds (see B96): 



Cpq Clp C\q 



(37) 



/C12 = /C03 — ► /C3 = 3^2CT 4 + ( — 2 + 18^2 —^2 — 4^3 



However, as can be seen from the above formulae, Eq.(p4)- 
(M), this is not true beyond the leading order (i.e, for the 
loop corrections in PT). 

In the limit ri2 — > 0, there is a simple correspondence 
between the 2-point and 1-point moments Idj — > K-i+j, or 
equivalently, cy — > Si+j. In particular, in this limit, one 
recovers known results for the 1-point moments in the SC 
model (see FG98, Appendix Al): 

AC11 = IC02 -* a% L = o 2 + (3 - 4// 2 + \v\ + «*) ° 4 + O (a 6 



+ 
and thus, 



V-2, + 6^2^3 + -^4 j CT 6 + O (c 8 ) 



C12 — > S3 = 3^2 + ( -2 + — v\ - 1v\ - 4^3 + -vA a 2 
C13 = C22 -> S4 = 12f 2 + 4^3 

+ (6 - 36^2 + 154 + 60i/f - 154 - 4^3 

— 20^2^3 — 6^2 V'i — 5^4 + 18^2^4 + 2^5) a , 

which provides a reassuring check for our results for the 2- 
point moments. 

3.2 Smoothing Effects 

Following FG98, the smoothing effects in the moments of 
the density field in cells of given resolution A, can be eas- 
ily incorporated. In particular, smoothing simply changes 
the v n coefficients of the SC model, Eq.(|12J), by introducing 
smoothing corrections in terms of derivatives of the logarith- 
mic slope of the rms fluctuation, 7 = dloga 2 /dlog A, 



1' 2 



I'-, 



V2 + ■ 



7 3 7 
V-& 1 V 2 7 H 

2 2 4 



4 28 

^4 — » C4 + -7-4i'27 + 2i'27+-V37- 



4 2 

o 7 



+ 



8 2 " 3 

r 27+ 27 7 - 



(38) 



where we have assumed for simplicity that higher-order 
derivatives of a are negligible (see also Juszkiewicz et al. 
1993). For an arbitrary power-law P(k), the above results 
can be trivially generalized (see Appendix A of FG98). 
Eqs.(|34|)-(|36[) along with the expression for the smoothed 
vertices v n , Eq.(^8[), can be directly compared to Nbody 
measurements (see §|4l below for a detailed discussion) . 

Alternatively, the PT solutions by B96 for the 2- 
point cumulants at tree-level (neglecting again higher-order 
derivatives of ex) read: n 



PT 
Cl2 



PT 
Cl3 



b8 7 
21 3 

11710 61 2 2 

1 7 H — 7 

441 7 ' 3 ' 



(39) 



Substituting the smoothed vertices, Eq.(p8h, into 
Eqs.(|34|)-(|36[), one sees that both results differ in the smooth- 
ing effects. For example, according to the SC model predic- 
tion, one finds 



(■12 



CIS 



68 2 

21 + 3° 


I 








11710 
441 


515 
42 


7 + 


17 
12 


2 
7 



(40) 



which shows that the 2-point skewness, C12, has an addi- 
tional smoothing correction 7/3 with respect to the PT re- 



S We have explicitly checked that the last term in Eq. (20) of B96 
(neglected in the above expression) makes little difference for the 
large values of ri2 where this expression applies. 
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suits derived by B96, Eq.fejf). Similarly, the 2-point kurto- 
sis, C13 , exhibits different smoothing effects in the SC model 
as compared to the PT results. Below we shall use Nbody 
simulations to investigate this discrepancy and discuss its 
interpretation. 



4 COMPARISON WITH SIMULATIONS 

In our numerical comparisons, we focus on two classes of 
models that approximate well basic observations of galaxy 
clustering and, in particular, the APM galaxy Survey (Mad- 
dox et al. 1990). 

The first class of models are Cold Dark Matter (CDM) 
ones. In particular, we use the flat lambda CDM model 
(ACDM hereafter) with a shape parameter V — 0.2 and 
£l m — 0.2 (f?A = 0.8). We generally use the outputs 
which have an amplitude of mass fluctuations, ug = 1, al- 
though we have also considered other outputs in the range 
(78 = 0.4 — 1.0. For completeness we also consider the "old 
standard" CDM (SCDM from now on) variant with F = 0.5 
and fl m = 1 (Qa — 0). Volume effects are assessed by 
comparing outputs from different box sizes in the range 
300 -600 ft" 1 Mpc. 

The second class of models considered have an APM- 
like linear power spectrum, so that, after evolution, it re- 
sembles closely the power spectrum inferred from the APM 
galaxy catalogue (see Baugh & Gaztanaga 1996). This model 
is normalized to <rg ~ 0.8 corresponding to the mean redshift 
in the APM catalog (z ~ 0.15). 

The small CDM simulations are the ones in Baugh, 
Gaztanaga & Efstathiou (1995) and the large ones (CDM 
and APM- like) are from Gaztanaga & Baugh (1998), where 
more details can be found. Errors are obtained from the 
dispersion in 10 (5) independent realizations of the small 
(large) simulations. An additional advantage of using such 
simulations is that their higher-order moments have been 
studied in detailed in the above cited references. 

Previous analyses of these simulations showed that the 
measured dark matter clustering, and in particular, its 1- 
point moments, are in good agreement with the leading order 
(also called tree-level) PT predictions, S q = £„/£ 2 ■ More 
recently, Fosalba & Gaztanaga (1998) found by using the 
SC model that such agreement extends further into non- 
linear scales when next-to-leading order PT predictions (the 
so-called loop corrections) are taken into account. 

We have also run a new set of simulations with CDM 
shape r = 0.25 but with Q. m = 1 (we call it Q m = 1 x V — 
0.25 CDM, where a CDM model with these values would cor- 
respond to a model with a low Hubble constant ft ~ 0.25. 
There is no need to make this correspondence, however: it 
can be considered to be an Q m = 1 model with a P(k) shape 
which is close to that seen in observations). These new sim- 
ulations can be used to test the sensitivity to fi m with in- 
dependent of the shape of P(k). They have a box size of 
300 ft" 1 Mpc and twice the particle resolution to the previ- 
ous sets (N — 200^ particles). Thus they can also be used to 
check for shot-noise and resolution effects. The computation 
of gravitational dynamics was started at z = 20 (instead of 
z — 10 in the other CDM sets) , so that they are less sensitive 
to possible Zeldovich Approximation transients (see Baugh, 




10 20 

r 12 (Mpc/h) 



Figure 1. The 2-point function £2 =< S1S2 > as a function of 
the separation r\% in Nbody simulations with a T = 0.2 CDM 
(triangles), T = 0.5 CDM (squares) and APM-like (circles) spec- 
trum. The continuous line shows the unsmoothed £2 predicted by 
linear theory in each model. The dashed line is the same linear 
prediction smoothed with the corresponding top-hat window. 

Gaztanaga & Efstathiou 1995, Scoccimarro 1998). There are 
5 realizations of the as = 0.5 output. 



4.1 Estimators 

The estimator of the density fluctuation at the ath pixel (or 
smoothed point) is 



5a 



Pa 

(p) 



1. 



(41) 



where p a is the mass density count in that pixel and (p) is the 
overall mean value. The estimator for the 2-point moments 
can then be generalized from Peebles & Groth (1975) into: 



=<S, 



>= 



1 



E 

a.b 



Stii W ab (r) 



(42) 



where N r = J^ o b W a b(r) is the number of pairs of pixels 
(or smoothed fields) at separation r in the sample, and the 
window function W a b(r) — 1 if pixels a and b are separated 
by \r a — ?b\ = rizdr, and otherwise. To obtain the reduced 
moments ktj we have to use the relations in equation (|22[). In 
paper II of this series we propose a variant of this estimator 
that produces better results for smaller samples. For the 
large periodic simulations we use here, the above estimator 
is simpler to code and works as well. 



4.2 The 2-point Function and the Linear Regime 

Figure [tj shows the comparison of the 2-point function 
^2(^12) with linear theory. The smoothing scale in this case 
are R f = 8 ft" 1 Mpc (LCDM, SCDM) and R f = 10 ft" 1 Mpc 
(APM-like). At large scales 7*12 > Rf, the 2-point function 
agrees well with the (smoothed) linear theory prediction. We 
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find that this is so for a wide ranges of smoothing scales, even 
when the smoothing radius is smaller (and the correspond- 
ing variance larger). This shows that the non-linear nature 
of small scale fluctuations does not affect significantly the 
large scale clustering, in agreement with PT predictions. 

On small scales, m ;$ Rf, the predictions follow closely 
the smoothed linear predictions (dashed line), but, as ex- 
pected, there are non-linear corrections when the smoothing 
scale is such that the variance is larger than unity. 

Note that the 2-point function in the APM-like model 
drops more sharply than the low-density model (see Fig |l|) . 
The SCDM model drop is even steeper. At large scales, the 
shape is similar to the APM-like but it crosses zero at a 
smaller scale. At the scale where the two point function be- 
comes zero the errors will be quite large, because of the 
relative effect of sampling fluctuations. This is also visible 
in the 2-point skewness, C12, which is subject to large biases 
at this point (see below). 

These models are quite different in shape and we will 
further explore the dependence in 7 by considering different 
smoothing scales. To see how the properties vary for a fixed 
scale and shape (7) we will also consider the statistics of 
different outputs (as parametrized by a&) from the same 
simulations. 

The agreement with linear PT shown in this section is 
a good test both for the code and estimator used. Higher- 
order 2-point moments are obtained from the same codes as 
the 2-point function, the only difference being the different 
powers of the density considered. 



4.3 The 2-point Skewness in the Weakly 
Non-linear Regime 

Figures 0-blshow the 2-point skewness, C12 = /Ci2/(£2(7" 2 ), in 
the weakly non-linear regime, that is when the smoothing 
scale is such that the variance is small: a ^ 1. The figures 
show a comparison with ACDM simulations for a wide range 
of values of a = 0.1 — 1.5 (a 2 = 0.02 — 2) corresponding to 
different scales, outputs and values of 7. The long-dashed 
line corresponds to the tree-level in PT prediction, Eq. (p9|) . 
The short-dashed line corresponds to the leading order pre- 
diction HI for the 2-point skewness in the SC: 



n , £2 

Cl2 = 2 V2 + v 2 — - 



21 3 r V21 3/ a' 2 



(43) 
Eq.©. 



while the dotted lines gives the limit £2 — * 0, e.g, 
The dot-dashed line includes the 1-loop SC correction (e.g, 
Eq. ( J34[ ) ) . Each panel is labeled with the smoothing scale Rf 
and the corresponding linear rms <?{Rf) and slope 7. 

In all cases we see there is a very good agreement with 
the leading order SC prediction. Note that in the limit 
£2 —>■ we recover the SC value in Eq. (|43i) (dotted line) 
rather than the rigorous PT prediction, Eq.(K9f) (long dashed 
line). This is a somewhat surprising result, as one would ex- 
pect the rigorous PT prediction to be more accurate than 
the SC approximation. A possible explanation could be that 
the regime of validity of PT is out of the dynamical range 
proved in our simulations. However this is unlikely to be the 



ACDM <7 8 =1.0 
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Rf= 11 a = 0.8 7=-1.63 
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Figure 2. The hierarchical coefficient C12 in the ACDM model 
as a function of the separation 7*12. Each panel shows results for 
different smoothing scale Rf as labeled in the Figure. Symbols 
show C12 in the erg = 1.0 output of the ACDM model, with errors 
corresponding to the variance in 5 realizations. The continuous 
line shows the leading order prediction for the skewness: S3 = 
Ci2(ri2 = 0) = 34/7 + 7. The long-dashed line corresponds to 
the PT result c 12 = 68/21 + 7/3 (B96V The short-dashed line 
corresponds to the SC prediction, Eq.rtla), which tends to C12 = 
68/21 + 27/3 (dotted-line) in the limit £2 -» 0. The dotted-dashed 
lines correspond to the 1-loop SC prediction, Eq.(B4). 
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■ Here by leading order in a PT calculation of the density fluc- 
tuation field, <5, we mean a term of order 6°, such as ^2/c 2 



Figure 3. Same as Fig.H for the trg = 0.5 ACDM output. 
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Figure 5. Same as Figured for the APM-like model. 



Figure 4. Same as Fig.0 for a 8 = 1.0 Cl m = 1 X T = 0.25. 



case as we have managed to measure C12 for scales as large as 
Rf <; 30 fo _1 Mpc and r\% <; 50 ft -1 Mpc. On those scales the 
variance is a 2 < 0.02 and £2 approaches zero, where pertur- 
bative expansions of 2-point moments should yield accurate 
predictions. 

On large scales the agreement with SC does not seem 
to depend on the value of the variance at the smoothing 
scale, except when we approach a ~ 1 and the perturbative 



r-0.5 CDM 



approach itself breaks down (see §4.4). 

Fig.W shows a comparison similar to Fig.H for the Q, m = 
1 x r = 0.25 model. The results are almost identical, indi- 
cating that they are insensitive to the value of Q m and to 
the particle resolution. 

Fig. Bl shows the corresponding comparison for the 
realistic APM-like model. It displays large fluctuations at 
ri2 — 50 ft -1 Mpc which reflects the fact that the 2-point 
function in the APM-like model goes to zero on that scale. 
A similar problem arises around ri2 — 40 ft -1 Mpc for the 
SCDM model in Fig n. As the two point function crosses 
zero quite sharply and the smoothing scale Rf is large, the 
zero crossing affects a large range of scales ri2 — 40 /i -1 Mpc 
(typically a range ±7?/). 

We can see Figured that the SC prediction gives a much 
better match to the simulations than the PT predictions, in 

Sreement with the cases discussed previously. Both in Fig 
ind |6| we can see how the SC model including the next-to- 
leading order (i.e, 1-loop correction) accurately reproduces 
the higher-order moments measured in simulations . 



4.4 The 2-point Skewness in the Non-linear 
Regime 

Figure M shows a comparison for smaller smoothing radius. 
Here the variance at the smoothing scale Rf ~ 3 — 4 is 
a 2 ~ 2. Notice that, as expected, significant departures 



3 
2.5 

2 
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g 
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Figure 6. Same as Figure [J for the T = 0.5 CDM model. 



from the tree-level prediction in the SC model are observed 
on such small scales, i.e, when the variance largely exceeds 
unity. This is in agreement to what is already known for the 
skewness 53. On the other hand, as can be seen in Fig m, the 
SC predictions work well on large scales, where £2 < 1 (e.g, 
for ri2 > 8/i _1 Mpc), even when a 2 > 1. This is similar to 
what we found for £2 in § 4.2 , which follows well linear theory 



even when the smoothing scale is fully non-linear. 

For small smoothing radius Rf the errors become much 
smaller and the SC gives a perfect match to Nbody mea- 
surements at large r"i2. This says that the SC reproduces 
well the non-linear transition between small and large scales 
found in numerical simulations. 

In general, for a ^ 0.8 the results agree quite well with 



© 0000 RAS, MNRAS 000, 000-000 



2-Point Moments I: Theory and Simulations 9 



AFM-like Rf=14.3 cj = 0.5 y = -l, 7 



3 4 



APM-like 



Rf = 3. (7=1.4 7 = -l.l 



ACDM 



.»■ ■-*. ■ai.&.W-W.^. 

Rf = 4.3 (7=1.5 7 = -1.3 



5 10 

r 12 (Mpc/h) 



50 



20 


1 ' 1 


i- 


15 
10 


- 


«\ 


--- A-~ 






- 


5 


I_ 




"" 


--■_. 


:-l-^"-" 


^ 




: 




1 




i 


20 


_ 


T : 


15 


- ~-_ 


-M 




^*~ 


i 


r 


10 






--:j|-- 


5 






i 





10 



20 40 

r 12 (Mpc/h) 



60 



Figure 7. Same as Figure 0, but for smaller smoothing scales 
(which means larger values of the variance). The top (bottom) 
panel shows the APM-like (r = 0.2 CDM) model. 



the perturbative SC model. In some cases, when 7 is large 
(flat slope) and a ~ 1 the amplitudes C12 tend to increase 
at all scales, including at large ria- This behavior is qual- 
itatively predicted by the perturbative SC expressions as 
given by Eq.(B4|), although in the non- linear regime, there 
is no reason to expect that the term ~ a 2 intrinsic to loop 
corrections in PT can accurately account for the gravita- 
tional evolution of the moments. We stress that the value of 
a, where deviations from SC prediction take place, depends 
on the value of the smoothing parameter, 7, in agreement 
with the general trend found in FG98 for the 1-point mo- 
ments: the larger the smoothing correction (or the smaller 
7), the smaller the non-linear contribution to the moments 
(see Table 2 in FG98). This can also be seen at large scales 
(<"i2 <; 8 ft -1 Mpc) by comparing the two panels in FigM 



Figure 8. The hierarchical coefficients C22 (top panel) and C13 
(bottom panel) in the APM-like model smoothed over a cell of 
Rf = 14.3 (<r = 0.5) as a function of the cell separation r\%. 
The continuous line shows the leading order prediction for the 
kurtosis: 54 = c\%(r\2 = 0) = C22(t*12 = 0). The long dashed 
line corresponds to the rigorous PT prediction. The short-dashed 
line corresponds to the SC leading order predictions, which tend 
to the dotted-lines in the limit £2 — > 0. The dotted-dashed line 
corresponds to the 1-loop SC prediction as given by Eqs.(p5|)-(pq). 



perturbative SC model, Eqs,(|35|)-(|36[), although noticeable 
differences appear when a ^ 0.5. This is due to the fact that 
the perturbative regime has a narrower domain of validity as 
one considers higher-order moments of the density field. This 
is so because loop corrections are relatively more important 
for the 2-point kurtosis (see the difference between short- 
dashed and dot-dashed lines in Figs H & H) than for the 
2-point skewness (see same lines in lower panels of Figs rJ & 



4.5 The 2-point Kurtosis 

The analysis presented for the 2-point skewness above, can 
be straightforwardly extended to higher-order moments of 
the density field. In particular, we shall compare the results 
in the SC model for the 2-point estimators of the volume av- 
eraged 4-point function, what we shall call the 2-point kur- 
tosis, Eqs.(B5j)-(pq), with numerical simulations. This will 
allow us to assess how general are the results found from 
the analysis of the skewness. Figure pi & B compare the SC 
prediction for the normalized 2-point kurtosis, C13 and C22, 
with the APM-like and CDM simulations, respectively. Note 
how the limit for large scales (dotted lines) is quite different 
for C13 — ► 6ff + 3f3 than for C22 — ► 4i/f. In both cases the 
SC predictions are in agreement with simulations. 

In general, for weakly non-linear scales, a ^ 0.5, the 
results agree quite well with the perturbative SC model. 
For larger values of a the amplitudes tend to increase at all 
scales, including at large ri2. 

The observed behavior is qualitatively predicted by the 



5 CONCLUSIONS 

We have presented perturbative results for the 2-point mo- 
ments of the dark matter density field within the SC model. 
In particular, we have derived expressions for the 2-point 
skewness and kurtosis in the weakly non-linear regime (see 

We have found an excellent agreement between the SC 
model and CDM and APM-like simulations at all relevant 
scales. When the variance a 2 5s 0.5 the perturbative SC 
model is in full agreement (within the error-bars) with mea- 
surements of the 2-point skewness and kurtosis from nu- 
merical simulations. When normalized to the second order 
statistics, we have confirmed that our results are insensitive 
to the cosmological parameters or to the amplitude of the 
initial fluctuations. We have tested the robustness of our 
findings against particle resolution and volume effects. 

However, we have observed that the domain of validity 
of the perturbative SC model gets narrower as one considers 
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Figure 9. Same as Fig. R for T = 0.5 CDM. 



as those we have explored here, so that in practice we should 
use these predictions (rather than B96) to compare to ob- 
servations (see paper II). Note that the proposal by Szapudi 
(1998) to determine the bias with cumulant correlators will 
not work in the regime we are exploring here (unlike B96 
predictions, the 1-point and 2-point SC predictions have the 
same smoothing corrections). 

The 2-point moments should provide a convenient tool 
to study the statistical properties of gravitational clustering 
for fairly non-linear scales and complicated survey geome- 
tries, as those probing the clustering of the Lya forest. In 
this context, the perturbative SC predictions presented here 
provide a simple and novel way to test the gravitational in- 
stability paradigm over a wide dynamical range. 
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higher-order moments of the density field, as seem in the 
analysis of the 2-point kurtosis (as compared to that of the 
skewness) . 

For larger values of the variance of the matter density 
field, a 2 , the perturbative approach is observed to break 
down, as expected. This is more extreme for smaller smooth- 
ing corrections. As we approach the strongly non-linear 
regime, a 2 <; 2, the leading order SC predictions match the 
simulations quite well on scales ri2 > Rf, while for ri2 ^ Rf 
there are important non-linear effects which cannot not ac- 
counted for within our perturbative approach, as one would 
expect. 

Our numerical simulations fail to match the leading or- 
der analytic predictions from perturbation theory (PT), de- 
rived by Bernardeau (1996, B96). In the unsmoothed case, 
these calculations are identical to the SC results presented 
here. But B96 found a different result for smoothed fields (ie 
compare Eq.[B9j with Eq.[M). Our failure to reconcile the 
simulations with the B96 predictions is puzzling. We have 
managed to extend our analysis to quite small values of the 
variance, i.e, a ^ 0.02, and we have considered large simula- 
tion boxes (up to 600 hT 1 Mpc) and also different resolution 
cells. It is unlikely that larger (realistic) simulations could 
resolve this problem as we are already probing scales where 
the two point function goes to zero and starts oscilating. On 
these larger scales, both simulations and predictions become 
uncertain. The deviations we find are quite significant, given 
the errors. One possible future approach could be to reduce 
the errors so as to be able to compare the predictions to 
earlier simulation outputs. This would allow us to explore 
the regime where a 2 ^ 0.02, while £2 is still positive. 

The idea that PT could fail for a 2 ~ 0.02 and start 
working for a 2 & 0.02 is possible, but unlikely. Experience 
with other PT calculations and comparison with simula- 
rtions (e.g., with the bispectrum, see Scocimarro etal 1998) 
indicates good convergence even for relatively large variance 
a 2 <; 0.1. It is reassuring that the perturbative SC approach 
works so well for a wide range of (realistic) situations, such 
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